      SUBROUTINE IMMDFL(A,NA,B,NB,H,NH,AM,NAM,BM,NBM,Q,NQ,R,NR,F,NF,P,N
     1P,IDENT,DISC,NEWT,STABLE,FNULL,ALPHA,IOP,DUMMY)
C 
C   PURPOSE:
C      Solve either the continuous or discrete time-invariant asymptotic
C      implicit (model-in-the-performance-index) model-following
C      problem.
C 
C   Subroutines employed by IMMDFL: ADD, ASYREG, EQUATE, LNCNT, MULT,
C      PREFIL, PRNT, SCALE, SUBT, TRANP
C   Subroutines employing IMMDFL: None
C 
      IMPLICIT REAL*8 (A-H,O-Z)
      DIMENSION A(1),B(1),H(1),AM(1),BM(1),Q(1),R(1),F(1),P(1),DUMMY(1)
      DIMENSION NA(2),NB(2),NH(2),NAM(2),NBM(2),NQ(2),NR(2),NF(2),NP(2),
     1IOP(1),IOPT(5),NDUM1(2)
      LOGICAL IDENT,DISC,NEWT,STABLE,FNULL,HIDENT
C 
      IF( IOP(1) .EQ. 0 ) GO TO 200
      CALL LNCNT(6)
      IF( DISC ) WRITE(6,25)
      IF( .NOT. DISC ) WRITE(6,50)
   25 FORMAT(/' PROGRAM TO SOLVE ASYMPTOTIC DISCRETE IMPLICIT MODEL-FOLL
     1OWING PROBLEM'//' PLANT DYNAMICS '/)
   50 FORMAT(/' PROGRAM TO SOLVE ASYMPTOTIC CONTINUOUS IMPLICIT MODEL-FO
     1LLOWING PROBLEM'//' PLANT DYNAMICS'/)
      CALL PRNT(A,NA,' A  ',1)
      CALL PRNT(B,NB,' B  ',1)
      IF( IDENT ) GO TO 75
      CALL PRNT(H,NH,' H  ',1)
      GO TO 100
   75 CONTINUE
      CALL LNCNT(3)
      WRITE(6,85)
   85 FORMAT(/' H IS AN IDENTITY MATRIX'/)
C 
  100 CONTINUE
      CALL LNCNT(4)
      WRITE(6,125)
  125 FORMAT(//' MODEL DYNAMICS'/)
      CALL PRNT(AM,NAM,' AM ',1)
      CALL PRNT(BM,NBM,' BM ',1)
      CALL LNCNT(4)
      WRITE(6,150)
  150 FORMAT(//' WEIGHTING MATRICES'/)
      CALL PRNT(Q,NQ,' Q  ',1)
      CALL PRNT(R,NR,' R  ',1)
C 
  200 CONTINUE
      N = NA(1)**2
      N1 = N + 1
      IF( .NOT. IDENT ) GO TO 300
      CALL SUBT(A,NA,AM,NAM,DUMMY,NH)
      CALL SUBT(B,NB,BM,NBM,DUMMY(N1),NB)
      GO TO 400
C 
  300 CONTINUE
      CALL MULT(H,NH,A,NA,DUMMY,NH)
      CALL MULT(AM,NAM,H,NH,DUMMY(N1),NH)
      CALL SUBT(DUMMY,NH,DUMMY(N1),NH,DUMMY,NH)
      CALL MULT(H,NH,B,NB,DUMMY(N1),NBM)
      CALL SUBT(DUMMY(N1),NBM,BM,NBM,DUMMY(N1),NBM)
C 
  400 CONTINUE
      IF( IOP(1) .EQ. 0 ) GO TO 500
      CALL LNCNT(3)
      WRITE(6,450)
  450 FORMAT(//' MATRIX HA - AMH')
      CALL PRNT(DUMMY,NH,0,3)
      CALL LNCNT(3)
      WRITE(6,475)
  475 FORMAT(//' MATRIX HB - BM')
      CALL PRNT(DUMMY(N1),NBM,0,3)
C 
  500 CONTINUE
      N2 = N1 + N
      N3 = N2 + N
      N4 = N3 + N
      CALL MULT(Q,NQ,DUMMY,NH,DUMMY(N2),NH)
      CALL MULT(Q,NQ,DUMMY(N1),NBM,DUMMY(N3),NBM)
      CALL TRANP(DUMMY,NH,DUMMY(N4),NDUM1)
      CALL MULT(DUMMY(N4),NDUM1,DUMMY(N2),NH,DUMMY,NA)
      CALL MULT(DUMMY(N4),NDUM1,DUMMY(N3),NBM,DUMMY(N2),NB)
      CALL TRANP(DUMMY(N1),NBM,DUMMY(N4),NDUM1)
      CALL SCALE(DUMMY(N2),NB,DUMMY(N1),NB,2.0D0)
      CALL MULT(DUMMY(N4),NDUM1,DUMMY(N3),NBM,DUMMY(N2),NR)
      CALL ADD(DUMMY(N2),NR,R,NR,DUMMY(N2),NR)
      IF( IOP(1) .EQ. 0 ) GO TO 600
      CALL LNCNT(3)
      WRITE(6,525)
  525 FORMAT(//' MATRIX ( HA - AMH TRANSPOSE)Q( HA - AMH)')
      CALL PRNT(DUMMY,NA,0,3)
      CALL LNCNT(3)
      WRITE(6,550)
  550 FORMAT(//' MATRIX 2( HA - AMH  TRANSPOSE)Q( HB - BM)')
      CALL PRNT(DUMMY(N1),NB,0,3)
      CALL LNCNT(3)
      WRITE(6,575)
  575 FORMAT(//' MATRIX ( HB - BM  TRANSPOSE)Q( HB - BM ) + R')
      CALL PRNT(DUMMY(N2),NR,0,3)
C 
  600 CONTINUE
      IOPT(1)= 0
      IOPT(2)= 1
      IOPT(3)= 1
      N5 = N4 + N
      CALL EQUATE(A,NA,DUMMY(N3),NA)
      CALL PREFIL(DUMMY(N3),NA,B,NB,DUMMY,NA,DUMMY(N1),NB,DUMMY(N2),NR,D
     1UMMY(N4),NF,IOPT,DUMMY(N5))
      IF(IOP(1) .EQ. 0 ) GO TO 700
      CALL LNCNT(3)
      WRITE(6,625)
  625 FORMAT(//' PREFILTER GAIN')
      CALL PRNT(DUMMY(N4),NF,0,3)
      CALL LNCNT(3)
      WRITE(6,650)
  650 FORMAT(//' MATRIX A - B(PREFILTER)')
      CALL PRNT(DUMMY(N3),NA,0,3)
      CALL LNCNT(3)
      WRITE(6,675)
  675 FORMAT(//' MODIFIED STATE VECTOR WEIGHTING MATRIX')
      CALL PRNT(DUMMY,NA,0,3)
C 
  700 CONTINUE
      CALL EQUATE(DUMMY(N4),NF,DUMMY(N1),NF)
C 
      IF( IOP(2) .EQ. -1000 ) RETURN
C 
      NF(1)=NB(2)
      NF(2)=NA(1)
      NP(1)=NA(1)
      NP(2)=NA(1)
      IOPT(1) = IOP(2)
      IOPT(2) = IOP(3)
      IOPT(3) = IOP(4)
      IOPT(4) = 0
      IOPT(5) = 0
      HIDENT = .TRUE.
      CALL ASYREG(DUMMY(N3),NA,B,NB,H,NH,DUMMY,NA,DUMMY(N2),NR,F,NF,P,N
     1P,HIDENT,DISC,NEWT,STABLE,FNULL,ALPHA,IOPT,DUMMY(N5))
      IF( IOP(1) .EQ. 0 ) GO TO 800
      CALL LNCNT(3)
      WRITE(6,725)
  725 FORMAT(//' GAIN FROM ASYREG')
      CALL PRNT(F,NF,0,3)
      CALL LNCNT(3)
      WRITE(6,750)
  750 FORMAT(//' SOLUTION OF ASSOCIATED STEADY-STATE RICCATI EQUATION')
      CALL PRNT(P,NP,0,3)
      CALL LNCNT(3)
      WRITE(6,775)
  775 FORMAT(//' EIGENVALUES OF P')
      NDUM1(1)= NA(1)
      NDUM1(2)= 1
      CALL PRNT(DUMMY(N5),NDUM1,0,3)
C 
  800 CONTINUE
      CALL ADD(F,NF,DUMMY(N1),NF,F,NF)
      IF( IOP(1) .EQ. 0 ) RETURN
      CALL LNCNT(4)
      WRITE(6,825)
  825 FORMAT(//' GAIN FOR MODEL-FOLLOWING CONTROL LAW, U = - F X  , F =
     1(PREFILTER) + (ASYREG)'/)
      CALL PRNT(F,NF,' F  ',1)
      N6 = N5 + NA(1)
      CALL PRNT(DUMMY(N6),NA,'A-BF',1)
      NDUM1(2) = 2
      N6 = N6 + N
      CALL LNCNT(3)
      WRITE(6,850)
  850 FORMAT(//' EIGENVALUES OF A-BF')
      CALL PRNT(DUMMY(N6),NDUM1,0,3)
C 
      RETURN
      END
